Skip to content

Note

Click here to download the full example code

topography from points

collar elevations sample the ground surface. topography triangulates them into a surface and a grid, clipped to the drilled area, and checks each collar against the surface through the others: a collar far off that surface was likely surveyed or typed wrong.

import sys
from pathlib import Path

HERE = Path(__file__).parent if "__file__" in globals() else Path.cwd()
sys.path.insert(0, str(HERE.parents[1]))
import boitata as bt
import matplotlib.pyplot as plt
import numpy as np
from common import HIGHLIGHT, map_axes, save

data = bt.datasets.stacked_sulphide_lenses()
collars = bt.PointSet.from_table(data["collars"], z="Z")
truth = data["topography"]

the delaunay surface is linear within each triangle and passes through every collar. by default the grid is clipped to a concave outline of the collars, buffered by their spacing; cells outside it stay null.

t = bt.topography(collars, cell=20.0)
z = t.grid["z"]
print(
    f"grid {t.grid.count[0]} x {t.grid.count[1]} cells of 20 m, {np.isfinite(z).mean():.0%} inside the outline"
)

true_z = truth.sample(t.grid.coords, "Z")
error = z - true_z
print(
    f"against the true surface: mean {np.nanmean(error):+.2f} m, 90% within {np.nanpercentile(np.abs(error), 90):.2f} m"
)

Out:

grid 63 x 70 cells of 20 m, 50% inside the outline
against the true surface: mean -0.18 m, 90% within 0.50 m

each collar's residual is its elevation minus the surface through the other collars. the default threshold, 3 robust standard deviations of the residuals, adapts to the data: on these clean collars it sits a quarter of a meter out and catches only the bends of the terrain between holes. a collar is flagged only where its residual also beats those of its neighbors, since one bad collar pulls its neighbors' residuals up too.

print(f"default threshold {t.max_residual:.2f} m: {t.flagged.sum()} of {len(collars)} flagged")
print(f"largest residual {np.abs(t.residuals).max():.2f} m")

Out:

default threshold 0.24 m: 12 of 289 flagged
largest residual 3.82 m

a typed elevation 10 m too high on one collar stands out against a threshold set from the survey's accuracy. the other collars past 2 m each have a neighbor with a larger residual, so only the mistyped one is flagged.

holes = data["collars"]["HOLE_ID"]
typo = collars.coords.copy()
typo[100, 2] += 10.0
checked = bt.topography(typo, cell=20.0, max_residual=2.0)
for hole, r in zip(holes[checked.flagged], checked.residuals[checked.flagged], strict=True):
    print(f"flagged {hole}: {r:+.2f} m")
print(f"other residuals above 2 m: {(np.abs(checked.residuals) > 2).sum() - checked.flagged.sum()}")

fig, ax = plt.subplots(figsize=(6.5, 6), layout="constrained")
bt.plot.section(t.grid, "z", axis="z", ax=ax, cmap="terrain")
ring = t.outline.parts[0]
ax.plot(*np.vstack([ring, ring[:1]])[:, :2].T, color="0.2", lw=0.8)
ax.scatter(*collars.coords[:, :2].T, s=4, color="0.2")
ax.scatter(*typo[checked.flagged, :2].T, s=60, facecolors="none", edgecolors=HIGHLIGHT, linewidths=1.5)
map_axes(ax, "Topography from collars, the mistyped collar circled")
save(fig, "topography")

Topography from collars, the mistyped collar circled

Out:

flagged DD0101: +9.80 m
other residuals above 2 m: 5

snapping collars to the surface

snap_to_surface moves each collar onto a surface, here the surveyed topography grid, and every hole moves with its collar. the report gives each hole's shift: the mistyped collar comes down 10 m, and the others move by the small gap between each collar and the gridded surface.

table = data["collars"]
collar = {c: table[c] for c in table.column_names}
collar["Z"] = typo[:, 2]
holes = bt.Drillholes(collar, data["surveys"])
snapped, report = bt.snap_to_surface(holes, truth, column="Z")
shift = report["shift"]
worst = np.argmax(np.abs(shift))
print(
    f"{report['hole'][worst]} moved {shift[worst]:+.2f} m; the other holes within {np.sort(np.abs(shift))[-2]:.2f} m"
)

Out:

DD0101 moved -10.02 m; the other holes within 0.09 m

Total running time of the script: ( 0 minutes 0.254 seconds)

Download Python source code: example_02_19.py

Download Jupyter notebook: example_02_19.ipynb

Gallery generated by mkdocs-gallery